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ABSTRACT 

Since globular clusters (GCs) are old, low- TV systems their dynamics is widely believed 
to be fully dominated by collisional two-body processes, and their surface brightness 
profiles are fit by King models. However, for many GCs, especially those with HST- 
resolved central regions, and 'extra-tidal' features, King models provide poor fits. We 
suggest that this is partly because collisionless dynamics is also important and con- 
tribute to shaping the cluster properties. We show using time-scale and length-scale 
arguments that except for the very centers of clusters, collisionless dynamics should 
be more imp ortant than collisional. We then fit 38 GCs analyzed by |Noyola fc GeT>] 
hardt ( 2006 1 with (collisional) King and (collisionless) DARKexp models over the full 
available radial range, and find that the latter provide a better fit to 29 GCs; for six of 
these the fit is at least ~ 5x better in term ofrms. DARKexp models are theoretically 
derived maximum entropy equilibrium states of self-gravitating collisionless systems 
and have already been shown to fit the results of dark matter N-body simulations. 
(We do not attempt fits with ad hoc fitting functions.) 
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1 INTRODUCTION 

Globular Clusters (GCs) are old and contain a relatively 
small number (~ 10 5 — 10 6 ) of stars, as compared to galaxies, 
and are therefore believed to have relaxed through two-body 
encounters, both strong and weak. In a Hubble time some 
GCs would have formed many compact binaries, reached 
core collapse, and possibly undergone several gravothermal 
oscillations ( jElson et al.|1987| [Meylan fc Heggie||1997 l. 

The theoretical treatment of this gravitational evolu- 
tion problem is not straightforward. While it is possible 
to describe realistic thermal equilibria for self-gravitating 
systems (see § [3| , the negative heat capacity of such sys- 
tems implies that these equilibria are not stable. Any re- 
duction in total energy leads to an increase in the kinetic 
temperature of the system, however, it is reasonable to as- 
sume that if GCs evolve slowly enough they go from one 
quasi-equilibrium state to another. These quasi-equilibrium 
states are often described by King equilibrium models, which 
attempt to take into account two key characteristics of 
GCs: relaxation through two-body encounters, and Galactic 
tidal field. To incorporate these characteristics, [King ( |1966[ ) 
starts with the Maxwell-Boltzmann distribution function, 
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f(E) oc exp([$o — E]/a 2 ), which corresponds to the isother- 
mal sphere of central potential $0 and constant velocity dis- 
persion a, and truncates it to imitate the result of spatial 
tidal stripping experienced by systems in th e Galactic tida l 
field, f(E) oc exp([$ - E]/a 2 ) - 1; see also |Michie| (| 1963 1. 



King-Michie models extend isotropic King models to include 
velocity dispersion anisotropy. 

It is generally believed that GCs are well-fit by King, 
or King-Michie models I Elson et al.|198 7 Meylan & Heggie 



19971. However, deviations from King models have always 



been known to exist, both at small and large radii. At small 
radii, the central density profiles, especially in HST data, 
sometimes lack the flat density cores required by King mod- 
els but instead exhibit a range of power-law cusps ( |Noyola| 



& Gebhardt 20061. At large radii, many clusters are bet- 



ter fit with power-law density profiles that remain shallower 
than King to the last observable poi nt (|McLaughlin &: van] 
|der Marel||2005| |Carballo-Bello et al.||2012[ ). Other GCs ap- 
pear to contain 'extra-tidal' stars, which make the surface 
brightness profiles more extended than King models. 

Because some GCs are known to show mass segregation, 
many authors argue that GCs should be fit with multiple 
King models, each representing some narrow stellar mass 
range. Profiles consisting of several (and as many as 10) 



mass components have been presented (Da Costa & Free- 
man|197"6 1 . Several authors have also pointed out that mod- 



els other than King can provide more conforming fits, like 
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Wilson models (Wilson 1975 McLaughlin & van der Marel 
2005), or the power law models | Carballo-Bello et al.|2012 |. 
The main drawback of using the latter two models is that 
they are ad hoc fitting functions. 

We are motivated by the widespread deviations between 
observations and King models and the wide variety of pro- 
posed fixes to ask whether King is really the best basic model 
for GCs. More specifically, we ask if the two-body encounters 
that involve one or few stars are as dominant a feature as 
it is generally believed, or if collective collisionless effects on 
scales of hundreds or thousands of stars are also important. 

There is a possible indication of collisionless dynamics 
in at least one cluster. WLM-1, studied in detail by |Stephens 
et al.| ( |2006[ ), has an ellipticity of e = 0.17 which they con- 
clude is not the result of rotation or Galactic tides. Instead, 
the shape is most likely due to anisotropic velocity disper- 
sions, as occurs in many ellipticals. Two-body effects would 
have wiped out differences in velocity dispersion in different 
directions, which implies that two-body relaxation has not 
been very efficient in this GC. 

Violent relaxation, a collisionless process, is brought 
about by fluctuations in the large-scale, or global potential 
of the system. The difference from collisional relaxation is 
that particles — stars — exchange energy with the large-scale 
potential, and not directly among themselves. In the case 
of galaxies, the potential fluctuations are usually the con- 
sequence of collapse and subsequent oscillations that occur 
during formation. In the case of GCs, the time since ini- 
tial collapse is generally longer than the two-body relaxation 
time-scale, and so additional causes for global potential fluc- 
tuations might be necessary if collisionless relaxation is to 
be relevant for present-day GCs, such as GCs elliptical or- 
bits in the Galaxy, disk passage, Galactic tides, and core 
collapse and bounce. For example, from the proper motion 



measurements of NCG 6397 Kalirai et al. (20071 conclude 



that the cluster has made frequent passages through the 
Galactic disk. 

Mixing is another aspect of collisionless relaxation; it 
does not result in energy changes of particles, but is also 
thought to be important in the formation of galaxies ( Mer- 



ritt 



2005). It is difficult to assign a specific mechanism re- 



sponsiblity, because even though the importance of potential 
fluctuations and mixing as drivers of collisionless relaxation 
is well established, the full understanding of relaxation and 
its detailed mechanics is still lacking. So instead of inves- 
tigating the relaxation mechanisms, we take a different ap- 
proach. As we describe in Section [3.2| it is possible to derive 
the final state attained by the system after collisionless re- 
laxation is complete, with no regard for the dynamics that 
lead to it. The final relaxed state can be described by the 
mass density distribution, which can then be compared to 
that of GCs. 

Therefore, in the present paper we do not undertake 
any dynamical experiments. We fit GC surface brightness 
(SB) profiles from the literature, with two functions: King 
and DARKexp, which represent fully collisional and fully 
collisionless systems, respectively. Section[3]shows that each 
model has been derived from maximum entropy arguments 
implying that these are effective quasi-equilibrium models 
for the final state of the two types of systems. We do not fit 
ad hoc or loosely motivated functions. 



2 ASSESSING COLLISIONALITY OF 
GLOBULAR CLUSTERS 

It is well established that two-body, or collisional processes 
are responsible for many properties of GCs, like the presence 
and hardening of binaries, mass segregation, and possibly 
energy transport and core collapse ( |Meylan fc Hcggic 1997[ ) . 

Here we ask, how important are collisionless effects 
in Globular Clusters? One way to address the question is 
to estimate relevant time-scales. The two-body relaxation 
time-scale for GCs is shorter than the Hubble time. How- 
ever, one must also compare it to the collisionless relaxation 
time-scale which proceeds on dynamical, or approximately 
crossing time-scales ( Binney fe Tremaine|1987 Lynden-Bell 
1967). 



The two-body relaxation time-scale is of the order, 
t2bd = (r/v)N/(8\nN), where r is the size of the system, 
v is the characteristic velocity of the particles, and N is the 
number of particles. The dynamical time-scale is, td yn = 
2n(r/v c ), where v c is the circular velocity and is of the same 
order as v. The ratio of the two, tdynjtibd = 167r(ln N)/N is 
2.3 for N = 100 and 6 x 10~ 3 for N = 10 s . This difference 
suggests that collisionless relaxation will drive the overall 
dynamical evolution of a GC, but interior to radii where 
N ~ 100, the two-body processes become more important 
for the evolution. 

Another way to access the relevance of collisionless ef- 
fects as a function of distance r from cluster center is to ask 
whether the cluster's acceleration field is grainy or smooth 
on the scale of individual stars. For this, one needs to esti- 
mate two length-scales. The first, r eq is the distance from 
the star where the typical acceleration induced by that star 
equals the acceleration due to the global potential; 



GM(< r) 



Gm 



(1) 



Here, the total mass of the cluster is related to the mass of 
a single star, and the total number of stars, Mtot = mNtot- 
The other length-scale, r sep , is the typical separation be- 
tween stars at radius r, r aep = [(4-7r/3)p(r)/m] -1//3 . The ra- 
tio of the two is 

r sep _ / 3M tot y/3/ M(<r) y/a 1/6 

r eq ~ Wp(ry*) \ M tot ) tot ( ' 

If r 3ep <C r eq , as one would find near the center of a 
GC, then the acceleration field and the potential are grainy, 
so the situation is collisional. 

Note that because this calculation assesses the graini- 
ness of the acceleration field (at distance r), it includes all 
types of collisional encounters — the strong ones with nearby 
stars and the weak ones with more distant stars. 

If r sep 3> r eq , for example, in the outer regions of a GC, 
then the acceleration field and the potential are smooth, and 
the situation is collisionless. 

It is interesting to see how the two collisionality crite- 
ria, one based on time-scales and the other based on length- 
scales compare. In Figure[T]the top panels show the results of 
the time-scale argument, and the bottom panel, the length- 
scale argument. The top panel plots the fraction of the en- 
closed mass, or the fraction of the total number of particles 
as a function of radius. The horizontal axis is in units of 
r_2, the radius at which the logarithmic space density slope 
7 = — dlog(p)/dlog(r), is equal to 2. Above, it was deter- 
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mined that central regions where N < 100 will be largely 
collisional. The horizontal dashed lines indicate this level 
for representative Ntot values of 10 3 and 10 6 . The bottom 
panel plots r se p/r eq vs. radius. Using the length-scale ar- 



Williams 2010 Barnes fc Williams|2012 \ is indistinguishable 
from the King DF. Madsen ( 1996 1 shows that when com- 



gument above, it was determined that r s 



divides 



collisional from collisionless regimes. In the bottom panel, 
that division is represented by the dashed horizontal line at 
r S e P = 3r eq . This specific relationship between r sep and r eq 
guarantees agreement between the time-scale and length- 
scale viewpoints for Ntot = 10 3 and 10 6 systems, and for 
DARKexp and King density profiles. The vertical dashed 
lines highlight that in both the top and bottom panels the 
collisional/collisionless transition occurs at the same loca- 
tion in the system, for a given size and model. 

We note that instead of t2bd as an estimate of the two- 
body relaxation time-scale one could have chosen the Spitzer 



mean time-scale (Spitzer 19871 which is based on the diffu- 



sion rate of particles through phase-space due to encounters, 

3 

t sp = 0.34— — . (3) 

G 2 mp In A 

Here, a, m, p and In A are the velocity dispersion, average 
stellar mass, mass density, and the Coulomb logarithm, re- 
spectively. For the DARKexp and King density profiles of a 
range of shape parameters, t sp is larger than t^bd by at least a 
factor of 10 over the relevant radial ranges. Using the Spitzer 
time-scale would not have produced an agreement between 
time- and length-scales that we have in Figure [l] but since 
t sp is always longer than t2bd it makes two-body effects even 
less important than we have estimated. 

We conclude that both time- and length-scale consider- 
ations lead to the same assessment of the collisionality of a 
system or part thereof, if t?,bd is used. And, regardless of how 
two-body relaxation is quantified, both approaches suggest 
that in a typical cluster with N to t > 10 4 , most of the body 
of the cluster, with the possible exception of the very center, 
would evolve collisionlessly. 



3 MODELS 

3.1 King-Madsen models 

Even though King models are physically motivated and have 
a full dynamical description provided by the King distribu- 
tion function (DF), there did not exist a rigorous derivation 
of that DF until |Madsen| ( |1996[ ). Madsen, like King and oth- 
ers before him, assumes that repeated two-body encounters 
in GCs require that one start with the Maxwellian velocity 
distribution and hence Maxwell-Boltzmann statistics which 
deal with classical particles that do not obey the exclusion 
principle. 

Madsen finds that the equilibrium state can be obtained 
as the maximum entropy state. Unlike the standard deriva- 
tion, which assumes that the occupation numbers in all en- 
ergy states are always large and hence Stirling approxima- 
tion is valid, he argues that no approximation should be 



made. Following Simons ( 1994 1, the exact form of the result- 
ing DF is derived, f(E) = [Ag exp(— f3E)], where [.] means 
round down to the nearest integer, and A, g are constants, 
and P is the constant inverse temperature of the system. 
This DF has discrete 'steps' that are especially pronounced 



bined with the Poisson equation, his DF gives density pro- 
files that are very similar to King's, and Hj orth fc Williams] 
(20101); iBarnes & Williamsl (120121) show that the smooth 



version of Madsen DF results in density profiles that are 



virtually identical to the standard King (1966) profiles. 

Note that in the case of the King DF, ' — 1' was added 
by hand to imitate the action of the tidal field. In Madsen's 
derivation, the corresponding modification of the DF results 
from a proper treatment of low occupation numbers. 



3.2 DARKexp models 

It has long been realized that the smooth appearance of ellip- 
tical galaxies, which consist of tens of billions of stars, cannot 
be due to two-body relaxation processes, like the rare strong 
encounters between stars, or the numerous weak encounters. 
Some faster acting relaxation mechanism had to be at work 
to ensure that the bulk of the ellipticals relax in well un- 
der a Hubble time. This mechanism would involve a rapidly 
changing potential of the system, which would induce indi- 
vidual stars to change energy. Because the mechanism has 
to act fast, it was called 'violent' relaxation by |Lynden-Bell| 
( |1967[ ). It is an example of collisionless relaxation, because 
particles exchange energy with the global potential and not 
through strong or weak 'collisions' among themselves. 

If we care only about the final state of the system 
and not the full dynamical history, we can use the tools 
of statistical mechanics. This approach to collisionless self- 
gravitating systems was originally put forward by |Ogorod-| 
|nikov|(| 19571), and iLynden-Belll \ 19671) , who argued that the 



final steady-state would correspond to the most likely state. 
|Lynden-Bell| | |1967| ) incorporated the collisionless nature of 
the stellar flow through an exclusion principle in phase- 
space: the collisionless Boltzmann equation states that a 
collisionless fluid is incompressible, hence phase-space ele- 
ments cannot be superimposed. In the non-degenerate limit, 
Lynden-Bell's theory predicted that the maximum entropy 
states would be isothermal spheres, as was also found by 



Ogorodnikov (19571. However, this result is unsatisfactory. 



While the maximization procedure imposed constraints of 
finite mass and energy, isothermal spheres are infinite. Fur- 
thermore, elliptical galaxies look nothing like isothermal 
spheres. Several ways out of this problem were proposed in 
the following decades. (See the |Hjorth fc Williams| ( |2010[ ) for 
a further discussion of the problems with the Lynden-Bell's 
approach and proposed solutions.) 

|Hjorth fc Williams| ( pOlO] ) argue that to apply statistical 
mechanics to collisionless self-gravitating systems requires 
one to make two important modification to the |Lynden-Bell| 
( |1967[ ) approach. First, the proper state space for collision- 
less systems, like dark matter halos and stars in elliptical 
galaxies, is energy space and not the standard phase-space. 
The rationale being that in a collisionless system in equilib- 
rium the particles' energies are fixed, so using energy space 
automatically ensures collisionlessness. Second, the expres- 
sion for the possible number of states, or entropy, involves 
n!, where n is the occupation number in energy space. When 
extremizing entropy, it is customary to simplify In n\ using 
the Stirling approximation, which is valid only for large n. 



at low occupation numbers. Its smooth version (Hjorth & 



Hjorth & Williams ( 2010 1 replace the Stirling approximation 
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with an expression that treats all n regimes very accurately; 
the implication being that low n regime, where n is 0, 1, 2, ... 
is relevant in self-gravitating systems. With these two mod- 
ifications they derive the most likely distribution in energy 
for the equilibrium systems, 

N(E) oc exp(-p[E - $ ]) - 1 = exp(0 o - e) - 1. (4) 

DARKexp is a single parameter family of models, with <f>o 
acting as a dimensionless potential depth. Density profiles 
based on DARKexp models can be found in Willi ams fc| 
|Hjorth| ( |2010| ). Unlike isothermal spheres, DARKexp sys- 
tems have finite mass and energy. A comparison with the 
results of collisionlcss dark matter simulations are presented 
in Williams et al. (20101. The simulations are a very good 
match to DARKexp models with 4>o around 4-5. 

We emphasize that both collisional and collisionless sys- 
tems, as described by the Madsen and DARKexp models, 
can be derived as the most likely statistical states. Aside 
from this, Madsen and DARKexp models have one other 
important feature in common: their derivation requires that 
the low occupation number regime is treated properly. The 
reason why this is the case for self-gravitating systems, while 
the Stirling approximation is adequate for other physical 
systems, is not yet definitively decided. What we can say is 
that realistic gravitationally bound systems have DFs that 
are truncated, leading to regions of state space that are very 
sparsely populated. It is therefore not surprising that en- 
tropy calculations, which depend on counting available en- 
ergy states, require those sparsely populated regions to be 
accurately accounted for. 



3.3 Our goals 

The primary goal of this paper is to find out whether GCs 
are better fit by King or DARKexp density profiles. These 
are both "first principles" models and not empirical fitting 
functions like the Wilson and power-law models. If a first- 
principles model is shown to fit the mass distribution of GCs, 
the most straightforward conclusion is that the physics that 
went into making the model applies to GCs. However, it is 
also possible that the goodness of fits is purely fortuitous, i.e. 
a complex combination of diverse dynamical effects happen 
to make the model a good fitting function. 

Our secondary goal concerns only those clusters that 
are well fit by the King model. Recall that the King and 
Madsen DFs differ somewhat because one is continuous and 
the other is discrete. It is not apriori obvious which one is 
appropriate for physical systems. The resulting King and 
Madsen density profiles are very similar, but not identical, 
so it makes sense to ask which of the two is preferred by real 
self-gravitating collisional systems. 



4 DATA 



Noyola & Gebhardt ( 2006 \ present the largest homogeneous 



set of non-parametrically estimated surface brightness (SB) 
profiles of 38 Galactic globular clusters. The inner regions of 
the profiles are derived from archival HST WFPC2 images, 
while the outer profiles are obtained from ground-based ob- 
servations. Estimating an unbiased smooth SB profile of the 
inner portions of GCs is not trivial, because a small number 



of very bright giant and horizontal branch stars introduces 
a considerable amount of shot noise. NG06 performed sev- 
eral simulations to determine the optimal way of overcoming 
shot noise, as well as photon noise. After extensive testing 
on synthetic data, they conclude that a certain combination 
of subtraction and masking of bright stars works best in re- 
covering their input profiles. The profiles of the outer radial 



regions were taken directly from Trager et al. (19951, who 



used Chebychev polynomial fits to the photometric points 
of ground-based data. The analysis of NG06 combines the 
HST and ground-based data into continuous composite pro- 
files of 38 clusters, which range from 2.36 to 3.68 decades in 
radius. For each cluster they publish 100 data points, spaced 
equally in log(r). 

The objectives of |Noyola fc Gebhardt] ( |2006 1 were to 
estimate the central SB slope of GCs, deproject the light of 
GC where possible, and estimate the central slope of the 3D 
light distribution. Hence, they derive uncertainties on the in- 
ner slope, but do not quote errors for the full radial range of 
the SB profile. It is probably reasonable to assume that the 
typical uncertainties in SB over the whole radial range are 
smaller than those in the individual HST data points, which 
are ~ ±0.15 mag/arcsec 2 . In this paper we calculate rms 
differences between NG06 fits and our theoretical models, 
and so do not use data uncertainties. Our rms values span 
the range from ~ 0.014 to ~ 0.8 mag/arcsec 2 (see Table [l]), 
which can be compared to the above quoted approximate 
uncertainty. 



5 FITTING GC SURFACE BRIGHTNESS 
WITH DARKEXP AND KING MODELS 

Dynamical studies of GCs are consistent with them having 
no dark matter. The main source of radial M/L variation 
in GCs is mass segregation. Even in the presence of mass 
segregation, radial M/L variations need not be large. Lane 
et al. (2009, 2010) find only small radial variations in M/L 
in their GC sample. In the radial range that the authors 
consider trustworthy their estimate of M/L varies by less 
than a factor of 1.5—2, and is consistent with being constant. 
Furthermore, the |Noyola &: G ebhardt ( 2006j) procedure of 
obtaining smooth SB profiles relies mostly on main sequence 
stars, and so minimizes the effects of mass segregation. In 
this work we assume that Mj L is constant with radius, and 
hence the shape of the SB profiles gives the shape of the 
radial mass distribution in these systems. One has to keep 
in mind that this is an approximation, and the true M/L 
must have some radial dependence. 

We use the non-parametrically smoothed profiles of 
NG06 as our input. We fit the SB of each Globular Clus- 
ter with 2D projected DARKexp and King models. DARK- 
exp models are characterized by a single shape parame- 
ter, (j>o, a dimensionless potential depth. King models are 
also characterized by a single shape parameter, <3?(0)/er 2 , 
a dimensionless combination of central potential and sys- 
tem's constant velocity dispersion. An alternate parame- 
terization is through the concentration parameter, ck = 
log(rtidai/r cor e); there is a monotonic relation between 
&(0)/a 2 and ck for King models. 

DARKexp models have a simple analytical expression 
for the energy distribution N(E), and King models have a 
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simple expression for the distribution function, f(E). How- 
ever, neither model has an analytical expression for the cor- 
responding density profiles, so these have to be obtained 
numerically, even for isotropic systems. We first calculate 
a library of DARKexp and King models and then compare 
them to the input profiles from NG06. 

Our fitting is done in the space of log(SB) vs. log(r), 
not the corresponding linear quantities, and the rms devia- 
tions between GCs and models are also calculated in the log 
space. While DARKexp and King models are one parameter 
families, we fit for three parameters, one shape parameter, 
and radial and SB normalizations. We use the full radial 
range for all 38 clusters presented in NG06, and give each 
of their fOO points equal weight. 

Table [I] presents our results. The core collapse status 
of the cluster is denoted by the labels c (collapsed) and c? 
(possibly collapsed), and were taken from NG06. In general, 
a cluster is classified as core collapsed if it has a cuspy, or 
steep SB profile extending to the center ( jMeylan fc Heggie] 
1997). In this paper we label clusters c and c?, but make 



no judgement as to their physical state. The second and 
third sets of columns, separated by vertical lines, show the 
DARKexp and King model fit parameters. DARKexp <j>o 
and King $(0)/ct 2 and ck are described above. In order to 
compare the two models directly, we have defined a new 
concentration parameter, cms and cki3 as C13 = log[r(7 = 
3)/r(7= 1)]. These concentration values are also shown in 
the Table. Next to the rms values for each of the two models, 
we quote, in parentheses, the separate rms for the inner and 
the outer regions of the radial profile. These regions each 
contain 50 of the total 100 points presented in NG06. These 
rms values are evidence that neither the inner nor the outer 
portion of the fit dominates the total rms. 

Figures [2] and [3] summarize our fits for DARKexp and 
King models, respectively. The horizontal axis is the poten- 
tial depth (shape parameter) of the model, and the vertical 
axis is the rms of the fit. The best-fit model for each of the 
38 clusters is represented by a solid dot, which is circled if 
the cluster is considered to be core collapsed. The best fits 
for individual clusters are shown in Figures |4|13| Red dashed 
curves are the NG06 data, blue curves are best-fit DARK- 
exp, and black are best-fit King profiles. The fit residuals 
(in A mag/arcsec 2 ) are shown in the bottom insets of each 
panel. Note that some residual curves show small amplitude 
fluctuations with radius, that look like noise. These come 
from the NG06 surface brightness data which is quoted to 
only 4 significant digits. This limited precision shows up as 
'noise', which is especially visible when the fits are good. We 
do not smooth over this noise here, but do smooth over it 
in Section [8] 

All results presented in this Section indicate that most 
clusters are better fit by DARKexp than King. Addition- 
ally, about a quarter of the clusters are fit very well, over 
the entire available radial range, by DARKexp. Figure [14] 
plots best-fit rms of DARKexp and King models; only eight 
clusters are better fit by King. Note that all eleven core col- 
lapse GCs are better fit by DARKexp than King models. 
Even if these are taken out, there is still a clear preference 
for DARKexp for non-core collapse GCs. Recall that the 
DARKexp model was built to reproduce the equilibrium 
states of collisionless systems, like dark matter halos. Ap- 
parently, many GCs are also well fit by it. We interpret this 



result to mean that the overall SB profiles of GCs may be 
less affected by collisional processes among its member stars 
than previously assumed. 

To further compare DARKexp and King fits, in Fig- 
ure [15] we show residuals from best-fits for all 38 clusters, 
for DARKexp (top panel) and King (bottom panel), respec- 
tively. The horizontal axis is normalized by the core radius 
of each cluster, which was defined by NG06 as the radius 
where the SB drops to half its central value. Note that this 
core radius definition is non-parametric and thus unrelated 
to the definition of the core radius in King models; the half- 
light radius was determined from non-parametric smoothed 
profiles. In addition to showing that DARKexp residuals 
are smaller, the figure also suggests that the King resid- 
uals show a systematic pattern, while DARKexp residuals 
are more random. If true, then this is an additional indica- 
tion that DARKexp are a better fit, at least typically, than 
King models. 

A skeptic might argue that the full radial range of a GC 
should not be fit with a single profile because different dy- 
namical processes are at play at small and large radii. The 
innermost radii might be affected by a central black hole, 
or core collapse and compact binaries, while the outermost 
radii may be populated by extra-tidal stars which are not 
in equilibrium with the global potential. We do not address 
these possibilities, aside from pointing out that these modi- 
fications would need to be described by multiple adjustable 
parameters. 

Before proceeding, we examine one more issue relevant 
for fitting, the radial extent of the data. The goodness of fit 
depends on the radial range being fit, and we have already 
mentioned that the shorter ranges available in the 1970's 
and 1980's were well modeled by the King profile. In Fig- 
ure |16[ the vertical axis shows the log of the ratio of the 
radial range given in NG06, and the horizontal axis the rms 
of the best-fit: solid dots represent DARKexp and empty 
triangles, King. The vertical solid and dashed lines in the 
top and bottom halves of the plot show the average rms 
for the DARKexp and King fits, respectively. As expected, 
a larger radial range leads to poorer fits for both models. 
DARKexp is moderately better than King for smaller radial 
ranges (bottom half of the plot). As the radial range gets 
larger (top half of the plot), DARKexp fits become notice- 
ably better fits than King. This too points towards DARK- 
exp being a better descriptor of GCs. 



6 KING VS. MADSEN MODELS 



Our secondary goal in this paper is to compare King (1966) 



and Madsen ( 1996[ ) models. The King and Madsen models 
are created by smooth and discrete distribution functions, 
respectively, and show corresponding differences in the den- 
sity profiles. Since it is not clear which of these two are more 
physically appropriate, we fit both models to twelve GC that 
are well fit, rms < 0.15, by King models (except NGC 6352, 
see the caption of Figure [3| . The clusters that are poorly fit 
by King tend to be much better fit by DARKexp (see Fig- 
ure 141, and so are unlikely to be better fit by the Madsen 



model. 

Figure [17] presents the results. The blue solid and 
dashed lines, and the solid and empty circles represent King 
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and Madsen fits, respectively. The best fits for the same clus- 
ter are connected by a red straight line segment. Eleven of 
the twelve clusters are better fit by the King model. Because 
the differences between King and Madsen profiles are small, 
we might have expected a breakdown closer to 50/50. The 
reason is probably the very sharp drop off in density that 
Madsen models have at large radii, a consequence of the 
discreteness of the distribution function at low occupation 
numbers. Real GCs do not show such steep drop off. 



7 DARKEXP VS. KING MODELS: A CLOSER 
LOOK AT SOME GCS 

In this Section we consider a few specific GCs and com- 
pare the standard description of their SB to the one using 
DARKexp. 

Modern data, especially HST and ground-based com- 
posite data sets such as the one provided byjNoyola fc Geb-| 
hardt (20061, have extended radial coverage and good res- 



olution down to very small radii. Here we illustrate, with 
one specific example, that a large radial range is often 
needed to discriminate between competing models. Consider 
NGC 6388, the cluster with the largest velocity dispersion, 
a 



18.9 km/s | |Pryor fc Meylan 1993) in the NG06 sampl 

~. — * o" +~ ncn" Itii; — Til: — ~,~-+ 



Using data spanning 3 to 250 , Illingworth & Illingworth 
(1976) fit a King model of c K = 1.75. Within these 2 radial 
decades, the King model fits the data well. If we truncate 
the NG06 data to the above radial range, our best-fit King 
model has ck = 1.77, consistent with earlier findings. How- 
ever, the full 3.14 radial decades clearly prefer DARKexp 
over King (Figure |9j, with rms values of 0.062 and 0.455, 
respectively. 

About 20% of Galactic GCs show deviations from King 
models by having steeper central SB profiles, i.e., cusps in- 
stead of cores. Moderately steep central cusps could be, in 
some cases, a signature of an intermediate mass black hole 
(IMBH), between 10 2 and 1O 4 M (though the presence of a 
shallow cusp is not a reliable indication of the IMBH; |Ves-| 
perini fc Trenti||2010 1. Steeper cusps are explained as post- 
core collapse clusters. However, as NG06 point out, clusters 
undergoing gravothermal oscillations should spend a small 
fraction of their lives in collapsed states. If all clusters that 
show cusps are assumed to be in a collapsed state, then we 
are catching a disproportionately large fraction of GC in this 
short-lived phase of their evolution. Consider the GC with 
the steepest cusp in the NG06 sample, NGC 6681 (M 70). 
The King fit has c K = 2.3 and rms = 0.194. The DARKexp 
fit has <f>o — 7.5 and rms — 0.065. The rather steep central 
cusp is very well accommodated by DARKexp. Whether it 
speaks for or against the cluster being core collapse or host- 
ing a central black hole is not clear, but the overall mass 
distribution over the entire 2.72 decades in radius is defi- 
nitely better accounted for by DARKexp. 



Another example is NGC 6715 (M 54). |Ibata et aL 
( |2009[ ) find density and velocity dispersion cusps within the 
central 0.3pc which they interpret as evidence of a 9400 Mq 
IMBH. However, they can also explain the cusp if the central 



stars have moderate radial anisotropy. Wrobel et al. (20111 



find no evidence for the IMBH in the X-rays, deriving an up- 
per limit on the Eddington ratio < 1.4 x 10~ 10 . Regardless 
of whether the cluster hosts an IMBH, DARKexp provides 



a better fit, including the central cusp, than King; the two 
rms are 0.104 and 0.595, respectively. 

Deviations from King models at large radii are also seen. 
These are generally attributed to extra-tidal stars, since GCs 
are expected to lose stars due to tidal stripping. Several GCs 
do show spectacular tidal tails. However, in some others, the 
'extra-tidal' distribution is circularly symmetric, and hence 
at odds with a tidally induced scenario. Such a distribution 



can be explained by evaporated stars (Kiipper et al. 2010 



20111. However, a DARKexp model also fits. For example, 



NGC 5694 was recently examined by Correnti et al. (20111 



using VIMOS / VLT. These authors note that they do not see 
a break, or any tidal tails beyond the cluster's tidal radius of 
3.15 arcmin. Instead, the stellar distribution smoothly con- 
tinues across the tidal radius, up to 10 arcminutes from the 
center. It shows an almost constant SB slope before and af- 
ter the tidal radius, which the authors fit with a power law, 
R~ 3 2 (see their Figure 3). In 3D this would be r -4 ' 2 , com- 
pletely consistent with DARKexp, which has an asymptotic 
outer slope of —4, independent of cf>o. Our best-fit DARK- 
exp model for NGC 5694 has <j> = 1.75, with rms = 0.076, 
compared to King model rms = 0.163. 

That truncation radii based on King model fitting may 
not correspond to tidal truncation due to the Galaxy is gen- 
erally acknowledged. In Figures [4 |l3| the empty downward 
triangles indicate the King model truncation radii taken 



from Gnedin & Ostriker (19971. The filled triangles indi- 



cate the truncation radius estimated from the strength of 
the Galactic tidal field, r 3 = r 3 (M GC /[3M(< r)]), where r is 
the Galactocentric distance, Mac is the GCs mass obtained 
from its absolute luminosity in V and a constant mass-to- 
light ratio of M/L = 2.5 in solar units, and M(< r) is the 
mass of the Galaxy interior to the current (i.e. not peri- 
center) location of the cluster. (Note that some truncation 
radii lie outside the limit of the figures.) 

Some GCs that are not well fit with King models and 
show signs of mass segregation are fit with multi-mass King 



models; see for example, |Da Costa & Freeman ( 1976 1. These 



more complex models are especially needed for clusters that 
cover a larger dynamic range in SB, up to 5 decades. NGC 



5272 (M3) is one example. Da Costa & Freeman ( 1976) show 



that the inner radial range can be fit with c = 1.29 King, 
while the outer with c = 1.98. Their actual fitted model 
has 10 mass ranges, for the data that span 2.4 decades in 
radius. NG06 data spans 3.6 decades, and the DARKexp fit 
has rms — 0.211, which appears comparable to the residuals 



of the 10- mass King profile (see Fig 3 of Da Costa & Freeman 



An an illustration of how different estimates of global 
parameters can be depending on the model fitted, we quote 
effective radii and enclosed mass for two representative clus- 
ters, which also happen to be the first two in Table 1. In 
the case of NGC 104, DARKexp and King models fit about 
equally well, while NGC 1805 is much better fit with DARK- 
exp. First, we note that we use half-mass radii instead of 
core radii because the latter are not defined for DARKexp 
models. Further, because the outer-most projected density 
profile of DARKexp falls off as r -3 , the mass enclosed scales 



as M(< r) 



This relatively slow fall-off makes "to- 



tal" mass a radius-dependent quantity, and requires one to 
decide where to truncate the DARKexp to calculate the "to- 
tal" mass. We chose to use the last NG06 data point. With 
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these choices, the ratio of the half-mass radii of best-fit mod- 
els for NGC 104 is Rdarkcx P / ' Ruing = 0.95, and the ratio 
of total masses is Mdarkbx P / 'Mxing = 1.0. For NGC 1851, 
the corresponding values are both 0.77. 



8 FITTING 3D LIGHT DISTRIBUTION OF 

GCC WITH DARKEXP AND KING MODELS 



Gebhardt ( 2006 1 who present smooth non-parametric com- 



Of the 38 GC in the Noyola & Gebhardt (20061 set, 12 



show central dips in the smooth non-parametric fits, and 
cannot be deprojected to yield 3D density profiles. For the 
remaining 26 clusters, NG06 present their non-parametric 
deprojections. We take these 26 GCs, smooth them to elim- 
inate the fluctuations due to the finite precision of their 
published results (see the sixth paragraph of Section [5] 
above), and compute 3D density slope, 7. These are plot- 
ted in Figure (18] For reference, we also plot five DARKexp 
and four King models, with <f>o = 5.66,4.0,2.83,2.0, 1.0 and 
Ck = 0.84, 1.25, 1.83, 2.35, respectively. These particular val- 
ues for cj>o and ck were chosen arbitrarily. 

Figure [18] is another way of looking at the same GC 
data. The 26 clusters are grouped by their DARKexp rms 
values, with the top left panel showing the 7 GCs which are 
fit by DARKexp very well, while the bottom right panel 
shows 6 GCs where DARKexp is a relatively poorer fit. 
Note that the bulk of the clusters have shallow inner density 
slopes, though in most cases these are not as shallow as those 
of King cores. DARKexp, on the other hand, has a range of 
slopes and slope derivatives at r < r_2 that match GC well. 
At r > r-2 some clusters are consistent with DARKexp pro- 
files, but a few show non-monotonic 7 behavior reminiscent 
in character of King profiles. Overall, even though DARK- 
exp provide better fits, similarities with King density profiles 
are also seen at some radii. 



9 CONCLUSIONS 



King (1966) models, which are based on the assumption of 
two-body collisional dynamics ( jMa dscn 1996), do not pro- 
vide good fits to surface density profile of many Globular 
Clusters, especially those whose data span a large radial 
range. Deviations occur at small and large radii, which are 
usually explained by physics beyond the King model. Moti- 
vated by these widespread departures, we ask if collisionless 
effects due to collective behavior of cluster stars could be 
important. In Section [2] we present two independent order- 
of-magnitude arguments that answer the question in the af- 
firmative. 

We then use previously created Globular Cluster sur- 
face brightness profiles to further test the relevance of colli- 
sionless effects. Our goal is made easier because there exist 
first-principles analytical models of collisional and collision- 
less systems, the King-Madsen and the DARKexp models, 
respectively. Both were derived based on maximum entropy 
statistical arguments. Both are one parameter families, and 
can be readily converted to density profiles and hence com- 
pared to the surface brightness profile of GCs (assuming GCs 
contain no dark matter and have constant mass-to-light ra- 
tios). We take the SB profiles from the work of Noyola & 



posite profiles based on HST and ground-based data and fit 
these with King, Madsen, and DARKexp density profiles. 

Our main finding is that DARKexp models fit consid- 
erably better than either King or Madsen models, over the 
entire available radial range. While this is not proof that col- 
lisionless dynamics has shaped the mass distribution of GCs, 
it does suggest it as an interesting possibility. It is already 
well known that the consequences of collisional processes, 
like the presence and hardening of binaries, mass segrega- 
tion and evaporation of stars are important in observed and 
simulated clusters (Meylan & Hcggic 1997). In this paper 
we demonstrate that the overall density profile is well de- 
scribed by a collisionless prediction, the DARKexp family 
of models. Taken together, these observations suggest that 
collisional and collisionless processes co-exist in GCs, but 
are responsible for different sets of properties. An alterna- 
tive explanation is that a good DARKexp fit is purely for- 
tuitous in that it allows for a cuspy central profile and so 
phenomenologically accounts for the effects of M/L varia- 
tion with radius as well as the effects of core collapse or post 
core collapse, etc., but the physics is not related to collision- 
less relaxation; DARKexp just happens to be a better fitting 
function to a very complex system. 

Because King models are apparently not the best 'typ- 
ical' model for GC, we caution against using these to derive 
cluster global structural parameters, like characteristic radii 
and enclosed mass, especially if these are used in other types 
of analyses, for example, to look for correlations with metal- 
licity, colors, age, etc. 

A secondary question we have addressed concerns an in- 
triguing property of self-gravitating systems. The derivation 
of Madsen and DARKexp models is based on statistical me- 
chanics; basically counting the number of particles in various 
states. In classical applications, all states have large occu- 
pation numbers, so the Stirling approximation is commonly 
used. Madsen and DARKexp models depend critically on 
the accurate treatment of the low occupation number regime 
and so raise the question of how the DF should behave when 
n — 1,2,3.... Should it be step-like, i.e. discrete, or should 



it be smooth? Madsen (1996) and King (19661 models are 
the two corresponding versions. (DARKexp is smooth.) A 
comparison of Madsen and King models for a subset of GCs 
shows the latter fit better. 

In the future it will be interesting to extend our analysis 
to include the dynamical information on Globular Clusters. 
Velocity dispersion profiles, though of poorer quality and 
of shorter radial range than surface brightness profiles, are 
available for some clusters, and can be used in combination 
with the surface brightness profiles to assess the state of the 
GCs ( |Zocchi et aLp012| |. 

Since DARKexp has already been shown to fit the re- 
sults of collisionless dark matter JV-body simulations quite 
well (both in the density, p(r) and energy space, N(E); 
Williams et~aL|2010| >, it is interesting to compare the best- 
fitting shape parameter, (fro, to those obtained in this work 
for GCs. A r -body generated dark matter halos have a nar- 
row range of <j>o around 4-5. GCs (excluding those that are 
poorly fit by DARKexp) mostly have <j>o between 1 and 3, 
because at approriate radii these DARKexp profiles have 
shallow, or flat density slopes. A few GCs have (j>o values 
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between 3 and 8. It is unclear why these values of 
preferred in either one of the systems. 
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-2 2 

log(r/r_ 2 ) 



Figure 1. A comparison of two order of magnitude collisionallity tests: top panel is based on t^ yn and t^hd time-scales, bottom panel 
based on length-scales. The two estimates agree. See Section [2] for details. 
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<j) Q of DARKexp fit 



Figure 2. The curves represent DARKexp fits to 38 GCs. The best-fit DARKexp model for each cluster is marked by a solid dot. Eleven 
clusters have best-fit <f>o = 0.05; these are not to be trusted because 0.05 was the smallest <j>0 value considered because the outer SB 
profiles of most of these clusters are steeper than what DARKexp can have, so no DARKexp model can provide a good fit. The circled 
dots represent core collapse GC. 
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Figure 3. Similar to FigureJS] but for King fits. The cluster with best-fit (0) / sigma? beyond 11 (at the right edge of the plot) is NGC 
6352; Figure [9] It has an upturn in SB at large radii, which is difficult to fit even with very large concentration King models. 
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12 3 1 2 3 



log(r/arcsec) 



Figure 4. DARKexp (blue) and King (black) best-fits to four GCs. The red dashed curve is the Noyola & Gebhardt (2006) data, which 
consists of HST and ground-based observations, jointed and smoothed by those authors. The surface brightness (SB) is in mag/arcsec 2 
in the V band. The rms between the data and models are indicated in the plot, together with the values of the dimensionless potential 
depth for DARKexp model, </>0i an d the King concentration parameter, ck ■ The fit residuals are shown in the bottom insets. The vertical 
scale is the difference in SB, in mag/arcsec 2 ; the its span varies between panels. The empty downward triangles indicate the King model 
truncation radii taken from Gnedin & Ostriker(1997). The filled triangles indicate the truncation radius estimated from the strength of 
the Galactic tidal field; see Section [7] 
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Figure 5. Similar to Figure [1] 
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Figure 6. Similar to Figure [3] 
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Figure 7. Similar to Figure [1] 
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Figure 8. Similar to Figure [1] 
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Figure 9. Similar to Figure [3] 
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Figure 10. Similar to Figure [4] 
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Figure 11. Similar to Figurc|4] 
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Figure 12. Similar to Figure [1] 
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Figure 13. Similar to Figure [4] 
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Figure 14. Rms of best-fit King and DARKexp models for 38 GCs. Circled dots represent eleven GCs that are considered to be core- 
collapsed. Overall, DARKexp models fit GCs considerably better than King models; out of 38 GCs only eight are better fit by King 
compared to DARKexp. 
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Figure 15. Residuals from best-fit DARKexp and King models for 38 GCs. Radius is normalized by the core radius defined by NG06 
as the radius where the SB drops to half its central value. This radius is unrelated to King core radius. In addition to being larger, King 
residuals also show a more systematic pattern compared to DARKexp residuals. 
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Figure 16. Radial range, \og(r max j 'r m i n ) of NG06 data vs. rms of the best-fit for 38 GCs. DARKexp fits are filled dots, King fits are 
empty triangles. The horizontal line divides the sample in to two roughly equal sets. In each set the solid vertical line gives the average 
rms for the DARKexp fits, while the dashed line give rms for King fits. 
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Figure 17. Similar to Figure [3] but only for 12 GCs where King models are a good fit. The solid lines and dots represent fits with 
King models, while the dashed curves and empty circles are the Madsen model fits. Red line segments connect best-fit King and best-fit 
Madsen models for the same cluster. 
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Figure 18. 3D log-log density profile slopes for 26 GCs deprojected by Noyola & Gebhardt (2006) (blue dashed) compared to DARKexp 
(black) and King (magenta) models. All the curves are normalized to have slope 7 = 2 at r = r_2. The GCs are grouped into four panels 
based on how well they are fit by DARKexp; rms ranges are shown in each panel. 
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Table 1. DARKexp & King model best-fit parameters for 38 GCs 
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